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It is widely accepted that the hippocampal place cells’ spiking activity produces a cognitive map 
of space. However, many details of this representation’s physiological mechanism remain unknown. 
For example, it is believed that the place cells exhibiting frequent coactivity form functionally inter¬ 
connected groups—place cell assemblies—that drive readout neurons in the downstream networks. 
However, the sheer number of coactive combinations is extremely large, which implies that only a 
small fraction of them actually gives rise to cell assemblies. The physiological processes responsi¬ 
ble for selecting the winning combinations are highly complex and are usually modeled via detailed 
synaptic and structural plasticity mechanisms. Here we propose an alternative approach that allows 
modeling the cell assembly network directly, based on a small number of phenomenological selection 
rules. We then demonstrate that the selected population of place cell assemblies correctly encodes 
the topology of the environment in biologically plausible time, and may serve as a schematic model 
of the hippocampal network. 
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I. INTRODUCTION 

The mammalian hippocampus plays a major role in spatial learning by encoding a cognitive map of 
space—a key component of animals’ spatial memory and spatial awareness 01 El- A remarkable property 
of the hippocampal neurons—the place cells—is that they become active only in discrete spatial regions— 
their respective place fields |[3l (Fig. [^). A number of studies have demonstrated that place cell activity 
can represent the animal’s current location HO, its past navigational experience |[6lEl, and even its future 
planned routes |[8l|3. Numerical simulations suggest that a population of place cells can also encode a 
global spatial connectivity map of the entire environment ODHEI- Hence, it is believed that the large-scale 
hippocampal representation of space emerges from integrating the information provided by the individual 
place cells, although the details of this process remain poorly understood. 

Experimental studies point out that the hippocampal map is topological in nature, i.e., it is more similar 
to a subway map than to a topographical city map ifTSUT^ . In lITOl we proposed a computational approach 
for modeling its structure based on several remarkable parallels between the notions of hippocampal physi¬ 
ology and algebraic topology. For example, Cech’s theorem asserts that the topological structure of a space 
X can be deduced from the pattern of overlaps between regions that cover it (for details see lITTl and Methods 
in ITOl). The argument is based on building a special simplicial complex TV, each fi-dimensional simplex of 
which corresponds to a nonempty overlap ofn+l covering regions, and demonstrating that the topological 
signatures of TV and X are same ifTTIl . Since the place cells’ spiking activity induces a covering of the envi¬ 
ronment by the place fields, called a place field map (Fig.[^), Cech’s theorem suggests that the place cells’ 
coactivity (Fig.[^), which marks the overlaps of the place fields, may be used by the brain to represent the 
topology of the environment. In |[T0l - [T^ it was demonstrated that place cell coactivity can in fact be used to 
construct a “temporal” analogue of the nerve complex, T, the simplexes of which, cr = [ci, C 2 , Ck], cor¬ 
respond to the combinations of coactive place cells, ci, C 2 , Ck (Fig. [IP)- Using the methods of persistent 
homology ||T8l [191 it was shown that the topological structure of T captures the topological properties of 
the environment, if the range of place cell spiking rates and place field sizes happen to parallel biological 
values derived from animal experiments ifTOlfTTl . 

However, it remained unclear whether it is possible to physically implement this information in the 
(para)hippocampal network. On the one hand, electrophysiological studies suggest that place cells showing 
repetitive coactivity tend to form so-called cell assemblies—functionally interconnected neuronal groups 
that synaptically drive a readout neuron in the downstream networks |EQ]-E3 —which may be viewed as 
“physiological simplexes” implementing T. On the other hand, the place cell combinations of T are much 
too numerous to be implemented physiologically. In a small environment, c.a. 1 x 1 m, thousands of place 
cells are active and the activity of 50-300 of them is near maximal level at every given location l[2T1l . The 
number of combinations of hundreds of coactive cells in an ensemble of thousands is unrealistically large, 
comparable to ~ 10^^^. The number of cells in most parahippocampal regions, which may potentially 
serve as readout neurons, is similar to the number of place cells Il24ll . This implies that only a small fraction 
of coactive place cell groups may be equipped with readout neurons, i.e., that the cell assemblies may 
encode only a small part of the place cell coactivities—those which represent a certain “critical mass” of 
spatial connections. 

Physiologically, the place cell assemblies emerge from dynamically changing constellations of synaptic 
connections and are commonly studied in terms of the synaptic and structural plasticity mechanisms ll25l - 
Ell. For a better qualitative understanding of the qualitative properties of the cell assembly network and 
practical modeling of the hippocampus’ functions, we propose a biologically plausible empirical approach 
that allows selecting the most prominent combinations of coactive place cells directly and demonstrate that 
the resulting population of cell assemblies is sufficient for representing the topology of the environment. 
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FIG. 1: Place fields and place cells. (A) The blue, green and brown dots, corresponding to the spikes produced by 
three different place cells, form well-defined spatial clusters, which represent their respective place fields. Spikes are 
positioned in space according to the animal’s coordinates at the time of spiking. (B) A place field map produced by 
an ensemble of 300 place cells with mean peak firing rate / = 20 Hz and mean place field size 5 ' = 14 cm located in 
a 1 X 1 m environment. (C) A short time segment of the spike trains produced by three place cells. The periods of 
the cells’ coactivity, marked by dashed lines, indicate overlap of their respective place fields (panel A): cells ci and 
C 2 are coactive in the region 12, cells ci, C 2 and C 3 are co-active in the region 123. (D) A simplex (T 123 represents 
schematically the spatial connectivity encoded by the coactivity of cells ci, C 2 and C 3 . Its ID edges correspond to 
pairwise coactivity, e.g. cr 12 represents the coactivity of cells ci and C 2 . 


II. THE METHODS 

Mathematically, the task of identifying a subpopulation of coactive place cell combinations corresponds 
to selecting according to biologically motivated criteria a subcomplex To of the full coactivity complex T. 
The cell assemblies correspond to the maximal simplexes of To, (i.e., the ones that are not subsimplexes 
of any other simplex), in contrast with the maximal simplexes of the coactivity complex, T, which can 
represent any largest combinations of coactive cells. The “cell assembly complex,” To, should satisfy 
several general requirements: 

I. Effectiveness. In the readercentric approach fill , each cell assembly drives a coincidence detector 
readout neuron in the downstream brain regions. Since the number of the readout neurons is comparable to 
the number of place cells, the total number of the maximal simplexes in To, A^maxCTo)^ should be comparable 
to the number of its vertexes, NdTo), 


A^max(ro) ^ NciTo). 

However, the algorithm for selecting To should reduce only the number of coactive place cell combinations 
and not the place cells themselves, meaning that the number of vertexes in T and in To should not differ 
significantly. In mathematical literature, the number of ^-dimensional simplexes of a simplicial complex is 
usually denoted as fk, and the list / = (/o, /i, , fd) is referred to as the complex’s /-vector |[30l. However, 
since in neuroscience literature the letter / is often used to denote firing rates, we denote the number of 
^-dimensional simplexes by A/. Asa shorthand notation, we use A^max to denote the number of the maximal 
simplexes and Nc the number of 0-dimensional simplexes in a given complex. 

II. Parsimony. To avoid redundancy, only a few cell assemblies should be active at a given location. 
Conversely, the rat’s movements should not go unnoticed by the hippocampal network, i.e., the periods 
during which all place cell assemblies are inactive should be short. 

III. Contiguity. A transition of the spiking activity from one cell assembly ct/ to another cr^+i occurs 
when some cells in ct/ shut off and a new group of cells activates in (T/+i (see Suppl. Movies). The larger is 
the subassembly crij+i = cri n cr^+i that remains active during this transition (i.e., the more cells are shared 
by cTi and (T/+i) the more contiguous is the representation of the rat’s moves and hence of the space in which 
it moves. The overlap between a pair of consecutively active simplexes can be characterized by a contiguity 
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FIG. 2: Topological loops: each horizontal bar represents the timeline 
of a topological cycle in T (T): OZ) loops (connectivity components) and 
the \D loops. Most cycles last over a short time before disappearing. A 
few remaining, persistent loops express stable topological information 
that may correspond to physical obstacles in the rat’s environment. The 
time required for the correct number of cycles to appear is interpreted 
as the minimal time Tmin required for the rat to learn the environment. 
The environment used in these simulations (Figure IB) is topologically 
connected (h^ = 1), and has one central hole (hi = 1). Thus, the topo¬ 
logical barcode of this environment—the list of Betti numbers (Zpq, ^i, 
^ 2 v) —is (1,1, 0, 0,). The last spurious loop (blue \D loop) disappears 
at about =4.6 minutes, which is the learning time in this case. 



index 




dim(cr,- n cr,+ i) 
Vdim(cr/) dim(cr/+i) ’ 


which assumes the maximal value ^ = 1 for coinciding cell assemblies and ^ = 0 for disjoint ones. In 
constructing a cell assembly complex, we expect that the mean contiguity over the simplexes in To should 
not be lower than in T. 

IV. Completeness. The cell assembly complex To should capture the correct topological signatures of 
the environment, such as obstacles, holes, and boundaries. For example, the lowest dimensional OD and 
ID loops in To represent, respectively, the piecewise and the path connectivity of the environment, as they 
are captured by the place cell coactivity. This information should emerge from the “topological noise” in 
a biologically plausible time period, comparable to the time required to obtain this information via the full 
complex, T (see ifTOlfTTl and Fig.|^. 


Place cell spiking is modeled as a time-dependent Poisson process with spatially localized rate 

{r-rc9 

Ur)-fee , 

where r is a point in the environment, fc is the maximal firing rate of a place cell c, and Sc defines the 
size of the corresponding place field centered at l[3ll . In a familiar environment, the place fields are 
stable, that is, the parameters fc, Sc and remain constant |[32l[33l. In our simulations, all computations 
were performed for ten place cell ensembles, each containing 300 neurons with an ensemble mean maximal 
firing rate of 20 Hz and a mean place field size of 30 cm. The place field centers in each ensemble were 
randomly scattered across the environment and most quantities reported in the Results were averaged over 
ten place field configurations. 

Spatial map. We simulated the rat’s movements through a small (1 x 1 m) planar environment (Fig.[^), 
similar to the arenas used in typical electrophysiological experiments (see Methods in lITOll ) over T - 25 
minutes—the duration of a typical “running session.” The spatial occupancy rate of the rat’s trajectory (i.e., 
the histogram of times spent at a particular location) and the frequency of the place cells’ activity are shown 
on Fig.[^,B. The mean speed of the rat is 20 cm/sec, so that turning around the central obstacle takes about 
7 seconds. 

By analogy with the place fields, we designate the spatial domain where a combination of place cells 
comprising a simplex cr is active as its simplex field, So- (Fig. |^). If the simplex corresponds to a cell 
assembly, then may also be referred to as the cell assembly field. Similarly to the place fields and 
the place field map (Fig|^), the collection of all simplex fields forms a simplex field map and the cell 
assembly fields form a cell assembly map ((Fig.[^). These maps provide a better “geometric proxy” for the 
rat’s cognitive map because they illustrate both the activity and the coactivity of the individual place cells 
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FIG. 3: Spatial maps. (A) Occupancy of spatial locations in a 1 x 1 m environment—a 2D histogram of the time 
spent by the animal in different locations. (B) Frequency of place cells’ spiking: each dot marks the location of a 
place cell’s center and indicates the corresponding appearance rate according to the colorbar. Higher appearance 
rates appear in the domain where the spatial occupancy is higher. (C) Simplex field map. The place field map for the 
same place cell ensemble is shown in Fig. (D) Spatial distribution of the frequency of the maximal simplexes’ 
appearances. Notice that, since place cells with higher appearance rates tend to produce higher order cell assemblies, 
which, in turn, have lower appearance rates, the spatial distribution of rates on B. and D. are complementary. 


((Fig.[^-D). In the following, the structure of these maps will be used to discuss our selection algorithms. 
If the distinction between a cell assembly map and a simplex map is not essential, it will be referred to as a 
space map. 

Population activity. To define the population code |[34l of place cell combinations, we construct place 
cell activity vectors by binning spike trains into 1/4 sec long time bins (for a physiological justification of 
this value see cniia). If the time interval T splits into n such bins, then the activity vector of a cell c is 

mc{T) = 

where specifies how many spikes were fired by c in the time bin. The components of normalized 
by the total number of spikes, define spiking probabilities, pc,k - O- A stack of activity vectors 

forms an activity raster illustrated on Fig. 

Two cells. Cl and C 2 , are coactive over a certain time period T, if the dot product of their activity vectors 
does not vanish, 

mc,{T) ■ nic^iT) ^ 0 . 

The component-wise or Hadamard product of two activity vectors 

^ci,c2 = Omc2 = [mci,imc2,i,mci,2mc2,2,,^ci,.z^c2,^] 

defines the coactivity vector of cells ci and C 2 , which can also be viewed as the activity vector of the 
corresponding ID simplex (Ti 2 = [ci, C 2 ], = ^ci,c 2 * Similarly, the Hadamard product of k vectors, 

^o-n...k = ^Ci,C2,...,Ck = ^C1 O^C2 O 

defines the activity vector of the simplex cr\2,„k = [^1,^2,, Ck\. 

For each activity vector, Mo-, we also define its bit array mapping into a binary appearance vector, Oo-, 
which indicates during which time-bins the corresponding simplex cr has made its appearance, i.e., Oo-j = 1 
iff Mo-j > 0. The appearance rate, fo-(T), of a simplex cr over a time interval T, is defined as the Li norm 
of its appearance vector, averaged over that time interval, 

U{T) = {\IT)Lia^^i. 

These appearance vectors and appearance rates allow distinguishing the intrinsic physiological characteris¬ 
tics of place cells’ spiking, e.g., their maximal firing rate, from the frequency with which these cells activate 
due to the rat’s movements through their respective place fields. While the maximal firing rate of a typical 
place cell is about 15 Hz |[T1, the frequency of their activation is much lower. 


Appearance rate, Hz 
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FIG. 4: An activity raster of a population of 20 
place cells over 250 time bins. Each row defines the 
activity vector of the corresponding place cell. The 
color of the ticks indicates the number of spikes con¬ 
tained in the corresponding bin, according to the col- 
orbar on the right. At every time step, the nonempty 
bins in the vertical column define the list of cur¬ 
rently active cells, i.e., the active simplex cr^. During 
the time interval T, marked by the two vertical blue 
lines, cell Ci6 is coactive with cig but not coactive 
with Ci 4 . 


III. RESULTS 


The simulated ensembles of 300 place cells in the environment shown on Fig. produced a coactivity 
complex?^ with about TVmax = 1000 maximal simplexes. Despite the high dimensionality of these simplexes 
(up to D = 35, mean D = 17), the characteristic dimensionality of a facet shared by two consecutively active 
simplexes, ct/ and (T/+i, is relatively low, so that the mean contiguity of T' is ^ = 0.6. This implies that, 
geometrically, if the simplexes of T are viewed as multidimensional tetrahedrons, the selected complex, 
To{ 6 ), assumes a highly irregular shape (Fig.[ST]A). 

More importantly, nearly 100% of the maximal simplexes appeared only once during the entire 25 
minute period of navigation, i.e., a typical maximal simplex’s appearance rate is low, fa- ~ 10“^ Hz. How¬ 
ever, a typical vertex activated about 200 times or every seven seconds, suggesting that some of the lower 
dimensional subsimplexes may be better candidates for forming cell assemblies. Is it then possible to build 
a cell assembly complex To by discarding the high-dimensional maximal simplexes with low appearance 
rates and retaining their subsimplexes that appear more frequently? We tested this hypothesis by identifying 
the combinations cr whose coactivity exceeds a certain threshold fo- > 6 , and studied the properties of the 
resulting simplicial complex as a function of 6 (Fig.[^). 

First we observed that, as soon as the appearance threshold is introduced {6 > 10“^ Hz), the high dimen¬ 
sional simplexes start braking up, releasing large numbers of lower dimensional subsimplexes: the number 
of ^-dimensional subsimplexes in a fi-dimensional simplex grows as combinatorial coefficient e.g., 
for n = n and k -1, C^g 44,000. As a result, the complex To{0) rapidly inflates. As 6 increases further 
(6 > 0.04), the number of “passing” simplexes decreases, and To( 6 ) begins to shrink in all dimensions (i.e., 
Nd(0i) > Nd( 02 ) for 6 i < 62 , for all D, (Fig.|^). Despite this, the number of maximal simplexes remains 
high: A^max = 30 x at 0 = 0.04 Hz, A^max - 7 x Nc 6 = 0.07 Hz and A^max = 3 x A^^ for the highest tested 
threshold, 0 = 0.1 Hz (Fig.|^), while their characteristic dimensionality drops from D - 17 to D = 7 at 
6 = 0.04 Hz and to D = 3 at 0 = 0.1 Hz. The mean contiguity index for this range of thresholds remains 
close to ^ = 0.7, indicating that the degree of overlap between the selected combinations of place cells is 
higher than in the original coactivity complex. 

However, raising the passing threshold 6 quickly destroys the geometric integrity of the resulting com¬ 
plex’s spatial map. As shown on Fig.[^ for 6 = 0.05 Hz, only ~ 50% of the environment is covered by the 
remaining simplex fields, and for 6 = 0.07 Hz the simplex map barely retains its one-piece connectedness: 


in some cases the complex To splits in two (the corresponding Betti numbers, bo, are listed in Table [ST 
for an illustration see Fig. [S3]). For 6 = 0.1 Hz, the complex fragments into multiple components (mean 
bo ~ 7) that are riddled with holes: the Betti numbers bn>o indicate the presence of hundreds of stable loops 
in higher dimensions. Thus, even if the coactive place cell combinations selected at 6 > 0.05 Hz could be 
supplied with readout neurons and would form cell assemblies, the resulting cell assembly network would 
not encode the correct spatial connectivity. 

An additional problem is that reducing the order of the assemblies violates the “assembly code” for 
spatial locations: every time several subsimplexes cri are selected from a high-order maximal simplex cr. 
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FIG. 5: A direct selection of the simplexes by appearance rates. (A) In the original coactivity complex T(6 = 
0), the maximal simplexes (Tmax appear on average but once during the entire observation period, resulting in low 
appearance rates (fo- < 10“^ Hz, blue line). Imposing four different thresholds 6 (color coded) raises the appearance 
rates of the selected maximal simplexes almost uniformly in all dimensionalities. (B) Cumulative distribution of the 
number of maximal simplexes A^ax over the selected simplexes’ dimension. In the T (6 = 0) case Amax exceeds the 
number of vertexes (Nc = 300, black horizontal line) by almost an order of magnitude. Small threshold values result 
in an explosive increase of A^ax which then begins to decrease for 6 > 0.04 Hz, remaining significantly higher than 
Nc for all four tested values of 6 . (C) The histograms of the maximal simplexes’ dimensionalities fit with normal 
distribution. The high mean dimensionality (D = 17) observed in the 7~(6 = 0) case reduces to 5 = 2.2 for ^ = 0.1 
Hz. The width of the distributions is about 50% of D. (D) The histograms of the number of the coactive maximal 
simplexes, fit to an exponential distribution demonstrate that the typical number of coactive simplexes is large, yS >10. 
All values are averaged over ten place field maps generated by ten place cell ensembles with the same mean peak firing 
rate and mean place field size. 


several overlapping simplex fields So-^ are produced in place of a single As a result, the parsimony of 
the representation is compromised: a location that was previously represented by a single simplex becomes 
represented by a few of its subsimplexes (Fig. |S^-B). Fig. 1^3 shows a histogram of the numbers of 
simultaneously active maximal simplexes in Tq: although most of the time only a few maximal simplexes 
are active, a coactivity of many of them (n > 25) is not uncommon. Conversely, while most of the time—on 
average 84% for the selected place cell ensembles—at least one simplex is active, longer inactivity periods 
are observed as described by double exponential distributed with the rate jS ^ 3.5 sec (Fig.|S2]). 

Overall, since most of the To-requirements listed in the Methods fail, we are led to conclude that the 
most straightforward selection rule, based on selecting high appearance rates, does not produce the de¬ 
sired tradeoff between the order of the assemblies, the frequency of their appearances, and the quality of 
topological representation of the environment. This failure motivates the search for alternative methods. 

Method I. To produce a more detailed approach to selecting coactive cell combinations, we observe that 
place fields are typically convex planar regions, and hence the existence of higher order overlaps between 
them actually follows from the lower order overlaps. According to Helly’s theorem, a collection of n > 
D + \ convex D-dimensional regions in Euclidean space will necessarily have a nonempty common 
intersection, if the intersection of every set of D -i- 1 regions is nonempty (see |[37l[38l and Fig.jS^). From 
the perspective of Cech theory, this implies that if n convex regions which cover a D-dimensional space 
contribute all the combinatorially possible D-dimensional simplexes to the nerve complex, then they also 
provide all the higher (up ion-\) dimensional simplexes to it. In a planar (D = 2) environment, this implies 
that a set of four or more place fields has a common intersection, if any three of them overlap. Moreover, 
although mathematically it is possible that three place fields exhibit pairwise, but not triple overlap, the 
probability of such an occurrence is low (Fig. |S^). A direct computational verification shows that if a 
triple of place cells demonstrates pairwise coactivity, then, in over 90% of cases, it also correctly encodes a 
triple spatial overlap. In other words, a “clique” of pairwise coactivities indicates the overlaps of all higher 
orders, which implies that the spatial connectivity graph Gj^ whose vertexes correspond to the place fields 
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FIG. 6: Selecting maximal simplexes via the pairwise coactivity threshold (Method I). (A) The appearance rates 
of the maximal simplexes computed for four different pairwise appearance rate thresholds 6 decrease as a function of 
their dimensionality. The values at D = 1 correspond to the value of the threshold imposed on the links’ appearance 
rate. (B) Cumulative distribution of the numbers of maximal simplexes, Amax, over the selected simplexes’ dimension. 
The numbers of cells Nc for each threshold value are shown by horizontal lines. The tendency of the maximal 
simplexes to outnumber the vertexes Amax > characteristic for small values of 6 , is reversed around 6 = 0.07 
Hz, where Amax and Nc level out. (C) The histograms of the maximal simplexes’ dimensionalities fit with normal 
distribution. The mean dimensionalities are similar to the ones produced by the previous selection method. The width 
of the distributions is about 50% of D. (D) The histogram of the number of coactive maximal simplexes, fit to an 
exponential distribution, shows that the expected number of coactive simplexes (yS ~ 4) is significantly lower than in 
the previous selection method. The procedure of averaging over the place field maps is the same. 


and links represent pairwise overlaps, encodes most simplexes in the nerve complex TV. 

As a reminder, a clique of an undirected graph is a set of pairwise connected vertices. From the combi¬ 
natorial perspective, a clique and a simplex have the same defining property: any subset of a simplex is its 
subsimplex and any subset of a clique is its subclique; a maximal clique is the one that is contained in no 
other clique. Hence, each graph defines its own “clique complex,” the ^-dimensional simplexes of which 
corresponds to the graph’s cliques with ^ + 1 vertices ||39l. 

The observation that the nerve complex induced from the place field map can be approximated by the 
clique complex of the place field pairwise connectivity graph, suggests that the corresponding coactivity 
complex T can also be built based only on pairwise, rather than higher-order, coactivities. This approach is 
well justified physiologically, since pairwise coactivity detector pairs of synapses are commonly observed 
001 El. The rule for defining the temporal analogue of Gj ^—the relational graph Gr —is straightforward: 
a pair of vertexes is connected in Gr if the corresponding cells c/ and Cj are coactive. Thresholding pairwise 
coactivity rates according to the rule 


fl fci,Cj ^ ^ 

|o if/c,c,<^. 


( 1 ) 


allows constructing a family of relational graphs Grio) over the pairs of place cells with high coactivity. 
The higher the threshold is, the sparser its connectivity matrix Cij and the smaller the number of maximal 
cliques and hence of maximal simplexes in the corresponding clique complex. Since in the following the 
graph Gat will not be used we will suppress the subscript “7^” in the notation for G^. 

We studied the relational graphs G( 6 ) and their respective clique complexes To(G( 6 )) = To( 6 ) as a func¬ 
tion of 6 . First, we observed that the appearance rates of the maximal simplexes in To( 6 ) become sensitive 
to the simplexes’ dimensionality (Fig.[^), implying that this selection procedure in effect attributes differ¬ 
ent thresholds to simplexes of different dimensions by using only one free parameter 6 . Second, the size 
of To( 6 ) is not as large as before. As shown on Fig. even for a relatively low threshold 6 = 0.05 Hz, 
the number of maximal simplexes exceeds the number of cells only marginally. For higher thresholds, this 
number steadily decreases: Nd(0i) < Nd( 02 ) for 6 i > 62 > 0.04 and D > 3 , though in lower dimensions 
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FIG. 7: Figure 7. Correction algorithms. (A) A spatial projection of the 2D skeleton of T(0) shows gaps and holes 
that compromise, respectively, the piecewise and path connectivity of T(0). If the links across the gaps and holes of 
T(6) are restored, then its correct connectivity structure may be regained. (B) A ‘hole’ produced by five connected 
vertexes is closed by restoring some of the previously discarded crosslinks. (C) A projection of the relational graph G 
into the environment, shown in grey. The edges added across the gaps are shown in red and the edges added to fill the 
holes are shown in green. The vertexes that are left disconnected due to low appearance rates of the edges connecting 
them to other vertexes are shown by blue dots. (D) A spatial map of the resulting “patched” 2D skeleton of T(6). The 
parameter values are = 15 , m/^ = 10 , and the lowered threshold for reintroducing the missing links is ^ = 50 . 


(1 < D < 3) this number may increase. The characteristic contiguity ranges between ^ = 0.65 at 6 = 0.05 
Hz to ^ = 0.72 at 0 = 0.14 Hz, which is higher than the value produced by the direct simplex selection 
method. Geometrically, this implies that the collection of maximal simplexes selected by pairwise thresh¬ 
old selection is more aggregated than the collection produced via direct simplex selection, i.e., the resulting 
complex To(6) is geometrically more similar to a “simplicial quasimanifold” (see Fig.[ST^). However, the 
number of place cells Nc drops as a result of discarding too many links with low appearance rate: Nc = 290 
aiO - 0.05 Hz and Nc = 100 atO - 0.14 Hz. At 0 = 0.1 Hz number of cells levels out with the number of 
maximal simplexes, A^max ~ A^c = 260. 

As before, raising the coactivity threshold degrades the spatial map. Ai6 > 0.07 Hz the simplex fields 
no longer cover the environment and at 0 > 0.1 Hz the map fragments into pieces (Fig. |S5)). However the 
resulting complex exhibits a much more regular topological behavior: the correct signature (Z7o = 1 , = 1 , 

^2 = 0, Z 73 = 0, ...) in To(6) appears at 6 = 0.05 Hz. The higher order Betti numbers (bn> 2 ) remain trivial 
at still higher 0s (Table [S2)\), even though the connectedness and path connectivity of the environment (bi 
and bo) become nfisrepresented. 

This improvement of the behavior of To(0) suggests that, despite all the shortcomings, the link-selection 
strategy may lead to a successful model of the place cell assembly network. After all, it is not surprising 
that a single selection rule does not resolve all the aspects of the cell assembly formation. Yet if it captures 
the essence of the process, it should be possible to correct or to adjust its outcome. For example, one of 
the difficulties faced by the coactivity selection algorithm is that, for high 0 , To(6) may brake into several 
pieces. However, the gaps between them are small. Thus, if a few discarded edges of the relational graph 
that originally bridged these gaps are retained, then the connectedness of 7o(0) may be spared (Fig. [^). 
Similarly, a “hole” in the relational graph is a linear chain of edges, connected tail to tail, with no shortcuts. 
However, if the links with the lower appearance rate (/ > 0/^, 6 ^ < 0) that span across the hole exist at 0 = 0, 
then they also can be restored (Fig. Iz^)- This may remove the non-contractible chains of ID simplexes 
in To(0) that compromised its path connectivity (Fig. |^,D). Thus, we implemented the following two 
rectification algorithms: 

1. Filling gaps: find pairs of vertexes Va and separated in G(0) by more than Ug edges and then 
test whether these vertexes are connected directly by links (from G(0 = 0)) whose appearance rate 
exceeds a lower threshold 6g < 0. If such links exist, add them to G(0) (red lines on Fig.[^). 

2. Closing holes: A closed chain containing m/^ > 4 edges in G(0), with no shortcuts, is likely to 
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A. B. C. 



FIG. 8: Statistics of the vertex degrees in relaitonal graphs. (A) The histogram of the vertex degrees k in the 
neighbor-controlled relational graph G(no), computed for four different hq (Method II) and fitted to a power law 
distribution P(k) ~ k~^. The graph demonstrates that G(no) is a scale-free network. (B) The same distribution on the 
log-log scale and an independent linear fit of the powers y. The confidence intervals of the two fits, ranging between 
+0.15 and ±0.3, overlap for each case. (C) In the pairwise coactivity threshold (Method I), the histogram of the 
relational graph’s vertex degrees is fit by negative binomial distribution, suggesting that G(6) is similar to a random 
network. 

produce a hole in To(6). We identified such chains and restored the discarded cross-links whose 

appearance rate exceeds a lower threshold 6h < 0. 

Thus, both rectification algorithms depend on two parameters: the length of the involved chains {ug for 
gaps and nih for holes) and the value of the reduced threshold 6g and 6h. In our numerical experiments, 
we found that the optimal value for the thresholds is 6h - Og - 50, and the parameters range between 5 
and 10 {nih) and 10 and 15 {Ug). Typically, each rectification procedure is applied once or twice before 
the right signature of To{6) is achieved, and this without producing significant changes of the complex’s 
structure, such as altering the appearance of its simplexes or increasing its size A^ax- As illustrated in the 
Table [S^, the correct signature in the “repaired” complex is achieved for all cases at 0 = 0.07 Hz. In 
particular, at 0 = 0.07 Hz we obtain a simplicial complex To with the correct signature, having Nc - 260 
vertexes and about the same number of maximal simplexes, Amax ~ A^. These maximal simplexes appear 
on average at a rate of fo- > 0.07 Hz, at least during every other run of the rat around the environment, and 
have dimensionality D = 6. As a result, the requirements to To are met and the maximal simplexes of To 
may represent hippocampal place cell assemblies that together encode a map of the environment, and hence 
To itself can be viewed as the “cell assembly complex.” 

Method II. A common feature of the appearance-rate-based selection rules is that the resulting simplicial 
complex reflects biases of its spatial occupancy: higher dimensional maximal simplexes concentrate over 
the parts of the environment where the rat appears more frequently. In particular, the relational graph shows 
a higher concentration of edges over the eastern segment of the environment (Fig.[^) where the occupancy 
rate is highest (Fig. |^). On the one hand, this is natural since the frequency of the place cells’ spiking 
activity certainly does depend on the frequency of the rat’s visits to their respective place fields, which 
therefore affects the hippocampal network’s architecture |[28ll2^ . In fact, this argument is at the core of the 
classical “hippocampus as a cognitive graph” model l|42l|43l, which proposes that the architecture of the 
hippocampal network is an epiphenomenon of the place cell coactivity. On the other hand, the physiological 
processes that produce synaptic connections may be more autonomous. For example, the CAS region of 
the hippocampus is anatomically a recurrent network of place cells whose spiking activity and synaptic 
architecture are dominated by the network’s attractor dynamics Il44f[46ll . 

These considerations lead us to test an alternative method of constructing the relational graph based on 
selecting, for every cell, its closest neighbors as defined by the pairwise coactivity rate fci,cj- Note that 
the resulting number of connections may be different for different cells: a cell ci may be among the hq 
closest neighbors of a cell C 2 , and hence ci and C 2 become connected, but the set of hq closest neighbors of 
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A. 



B. 





FIG. 9: Figure 9. Selecting maximal simplexes via best neighbor selection (Method II). (A) The appearance rates 
of the maximal simplexes in the simplicial complex 7“(^o), computed for four different values of hq (color coded), 
decrease as a function of their dimensionality. (B). Cumulative distribution of the number of maximal simplexes A^ax 
over the selected simplexes’ dimension. For the tested values of hq, the fixed number of vertexes Nc = 300, indicated 
by the horizontal black line, is close to the number of maximal simplexes. For hq = 1, the values Amax and Nc come 
closest. (C). The histograms of the maximal simplexes’ dimensionalities, fit to the normal distribution, indicate that 
for the relational graph with a similar number of links, the mean dimensionalities of the maximal simplexes are smaller 
than in in the complex built via the threshold-selection method. The width of the distributions is about 40% of D. (D) 
The histogram of the number of coactive maximal simplexes, fit to an exponential distribution. An expected number 
of coactive simplexes ranges between jS = 2 and 13 = 6. The procedure of averaging over the place field maps is the 
same. 


a cell Cl may not include C 2 , which bears a certain resemblance to the preferential attachment models (471. 
As a result, the vertex degrees k of the (undirected) relational graph may differ from one another and from 
riQ. A direct computational verification shows that k is distributed according to a power law, P(k) ~ k~^, 
where y ranges, for different hq, between y ~ 2 and y ~ 4 (Fig. [^, which implies that G(no) demonstrates 
scale-free properties (4711481 characteristic of the hippocampal network (43|50l. In contrast, the histogram 
of the vertex degrees in the threshold-controlled relational graph G(6) may be fit with the negative binomial 
distribution (Fig.[^), which indicates that G(6) is similar to a random graph. 

This neighbor-selection method for building the relational graph G(no) has a number of other immediate 
advantages over the threshold-controlled construction of G(6). For example, no cells are excluded from To 
due to the low appearance of the edges connecting to them. As a result, the simplex fields are distributed 
more uniformly (Fig.|S6|), which helps capture the correct piecewise connectedness of the environment. 

By studying the properties of the clique complexes produced by the relational graphs G(no) for no = 2, 
4, 7 and 12—parameters chosen to produce similar numbers of edges as in the previous method—we found 
that the number of maximal cliques in To(G(no)) is typically lower than in To(G(6)). The appearance rates 
of maximal cliques in G(no) are more scattered and less sensitive to dimensionality than in G(6) (Fig. 
and Fig. [S6| ). The number of maximal simplexes in To(no) remains close to the number of cells (Fig. [^) 
and their dimensionality is lower than in the threshold-based selection approach (Fig. [^). The contiguity 
index in all complexes ranges between to ^ = 0.67 and ^ = 0.71. The coverage of the space with the 
simplex fields improves with growing no (see Fig.[S6|) —for no > 2 the complex To(no) is connected, while 
the behavior of bo is more regular (see Table [S3]). However, the path connectivity of the complex To(no) 
remains deficient for all no because the number of stable spurious ID loops remains high (Table [S3]). After 
filling the gaps and closing the holes, most complexes constructed for no > 1 acquire correct topological 
signatures (Table [S3]), and the requirements to To are satisfied. Thus, the simplicial complex obtained by 
the neighbor selection method for no>l can also be viewed as a “cell assembly complex,” meaning it can 
serve as a formal model of the place cell assembly network with mean contiguity ^ = 0.7. 
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IV. DISCUSSION 

The proposed approach allows creating schematic models of the place cell assembly network—the cell 
assembly complex To—by controlling the basic phenomenological parameters of place cell coactivity and 
then relating the network’s architecture to the net topological information it encodes. The coactivity com¬ 
plex T was previously used for representing the pool of place cell coactivities |[T0l - IT^ . Specifically, the 
low order (pair and triple) coactivity events were used to construct the 2D skeleton of T, and then its OD 
and ID topological loops were matched with the corresponding topological signatures of the environment. 
In the current study, the entire pool of place cell coactivities is used to model the full place cell assembly 
network, including the higher order assemblies representing both the low-dimensional spatial environment 
as well as high-dimensional memory space ll2Tl[5T1l . The learning times Tmin estimated from the dynamics 
of the OD and ID loops in To remain close to the learning times computed for the full coactivity complex T 
(see Table [S^. This implies that the selected, “core” pool of coactive place cell combinations captures the 
topological structure of the environment as fast and as reliably as the entire set of the place cell coactivities. 

We view the proposed algorithms as basic models of a more general “phenomenological” approach, 
one which can be further developed along several broad lines. First, the structure of the relational graph is 
currently deduced from the activity vectors defined over the entire navigation period T - 25 minutes. A 
biologically more plausible selection algorithm should be adaptive: the structure of the relational graph at 
a given moment of time t < T should be based only on the spiking information produced before t. Hence, 
in a more advanced model, the structure of the relational graph should develop in time, and in general 
the cell assemblies comprising To should be derived using synaptic and structural plasticity mechanisms. 
Second, the selection criteria in Methods I and II above may be individualized: the appearance threshold 
used to construct the relational graph can be assembly-specific, i.e. 6 - 6{cr), so that the properties of 
the resulting network would be described in terms of the probability distribution of the threshold values 
across the cell assembly population. Similarly, the number of closest neighbors can be made cell-specific, 
fio = noici), which should permit better control over the topological properties both of the network and of 
the cell assembly complex. Third, threshold control can be implemented using different coactivity metrics, 
for instance via the pairwise correlation coefficient 

p(ci,C 2 ) = (mi • m 2 )/|mi||m 2 |, (2) 

which would connect cells with correlated spiking (irrespective of their firing rates), in contrast with the 
metric ([^, which does the opposite. In general, two metrics p and p', produce relational graphs with differ¬ 
ent topologies. Nevertheless, they may produce similar or identical large-scale effects, such as generating 
topologically identical cell assembly complexes Tq, or exhibit similar learning times, Tmin. Identifying 
classes of metrics that produce topologically similar results will be examined in future research. 
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VII. SUPPLEMENTARY FIGURES 


A. B. 



FIG. SI: Simplicial complexes. (A) A schematic representation of an irregular simplicial complex in which 
the number of maximal simplexes is larger than the number of vertexes (black dots). The maximal \D simplexes 
are shown as blue segments, the 2D simplexes as gray triangles and the 3D simplexes as pink tetrahedrons. (B) A 
simplicial “quasi-manifold,” (3, which has a similar number of vertexes and maximal simplexes of different dimen¬ 
sionalities. If each maximal simplex, e.g. (1,2,3) or (10,7, 8,9), corresponds to an assembly of place cells driving a 
readout neuron (red dots), then (3 is a cell assembly complex. Dotted lines represent synaptic connections from the 
place cells to the readout neuron. 








16 


A. 



0 = 0.05 Hz 


2000 4000 6000 

Number of max simplex 



0 = 0.07 Hz 


0 = 0.1 Hz 


110 


5 o' 


2000 4000 6000 

Number of max simplex 


2000 4000 6000 

Number of max simplex 



112 


6 




FIG. S2: The simplicial complexes To(6) constructed by direct selection of coactive combinations, for four 
different values of 6 . (A) The appearance rates of simplexes, arranged from left to right according to their dimension. 
Each dot corresponds to a maximal simplex whose dimension is color-coded according to the colorbar on the right. 
(B) Spatial distribution of the dimensionalities of the selected simplexes. (C) Spatial distribution of the appearance 
rates of the selected simplexes. (D) The histograms of the lapse times, fit to double exponential distribution (blue 
line), and the value of the fitted distribution’s rate (3. (E) Spatial projections of the 2D skeletons of the Tq{6 ). Data for 
all panels is computed for a specific place field map for illustrative purposes. 
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D. 


bo= 1, b, = 0 


E. 



bo= 2, b,= 1 


F. 


bo= 2, b, = 2 




FIG. S3: Low-dimensional topological features of the cell assembly complexes captured by place field map. 

(A) A place filed map corresponding to a singly-connected (/^o = 1, marked by the dot) complex with no non- 
contractible loops {h\ = 0). (B) A place field map corresponding to a complex with correct list of Betti numbers 
(correct topological barcode (HI): the physical hole (see Figures 1 and 3) produces one non-contractible loop (red 
circle, hi = 1). (C) A map containing a spurious hole in the upper-right comer produces an extra persistent loop (the 
additional small circle, net = 2). (D) A map containing two disconnected pieces marked by the black dots (net 
Z^o = 2) and having no non-contractible loops {b\ =0). (E) A map containing two pieces and one persistent spurious 
\D loop. (F) The green and the brown place fields at the top connect, yielding another persistent spurious \D loop. 
Compare these illustrations and the topological barcodes to the illustrations and Suppl. Movies provided in (5^ 
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FIG. S4: Schematic illustrations of spatial maps. (A) Schematic representations of 3^ order simplex fields, a, b, 
c and d, encoded by third order maximal simplexes, = (vi, V2, ve), (Tc = (v2, V3, V5), etc. (B) If the 2D simplexes 
are discarded and their \D faces are retained, then the second-order simplex fields are produced, shown here as 
overlapping shaded regions. The original simplex fields a, b, c and d are now represented by the coactivity of three 
pairs, e.g., a is represented by = (vi, V 2 ), cra ,2 = (v 2 , ve) and cra ,3 = (vi, V6). (C). The three place fields on the left 
exhibit pairwise, but not triple overlap. In the generic spatial configuration shown on the right, pairwise overlapping 
place fields also produce a triple overlap. (D) Four pairwise overlapping convex regions in 2D produce all the higher 
order (triple and quadruple) overlaps. 
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FIG. S5: Simplicial complexes Tq{6 ) constructed via pairwise coactivity selection (Method I) for four different 
threshold values. (A) The appearance rates of the maximal simplexes, arranged according to their dimension, demon¬ 
strate remarkably tight, graph-like distribution. The color of the dots corresponds to the dimension of the simplexes, 
as indicated by the colorbar on the right. (B) Spatial distribution of the dimensionalities of the selected simplexes. 
(C) Spatial distribution of the appearance rates of the selected simplexes. (D) The histograms of the lapse times, fit 
to double exponential distribution, and the value of the fitted distribution’s parameter (3. (E) Spatial projections of the 
2D skeletons of Tq{6 ). Data for all panels is computed for a specific place field map for illustrative purposes. 
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FIG. S6: The selected simplicial complexes 7o(^o) constructed via the closed-neighbor selection algorithm 
(Method II) for four different values of hq. (A) The appearance rates of simplexes, arranged according to their 
dimensions, color-coded as indicated by the colorbar on the right. (B) Spatial distribution of the dimensionalities of 
the selected simplexes. (C) Spatial distribution of the appearance rates of the selected simplexes. (D) The histograms 
of the lapse times, fit to double exponential distribution, and the value of the fitted distribution’s parameter (3. (E) 
Spatial projections of the 2D skeletons of 7o(^o)- 
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VIII. SUPPLEMENTARY TABLES 


e = 0.04 

Hz 

i 

9 = 0.05 

Hz 

II 

p 

b 

Hz 

e = 0.1 Hz 

1 

139 

509 

569 

1 

171 

393 

410 

2 

198 

324 

214 

5 131 108 15 

1 

126 

409 

480 

1 

175 

353 

286 

1 

203 

247 

144 

7 147 71 17 

1 
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853 

2 

211 
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2 
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330 

419 

6 139 152 41 

1 
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1 
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369 

333 

1 

201 

296 

152 

7 139 50 7 

1 
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422 
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2 

172 

342 

453 

2 

185 

312 

320 

13 122 107 73 

1 

130 

485 

775 

1 

186 

447 

392 

1 

190 

335 

250 

4 164 84 27 

1 

106 

459 

655 

1 

202 

352 

603 

1 

177 

395 

353 

8 179 117 25 

1 

120 

519 

787 

1 

194 

397 

449 

1 

215 

333 

237 

8 160 109 63 

1 

131 

441 

757 

1 

165 

459 

499 

1 

203 

339 

197 

8 149 65 15 

1 

160 

364 

575 

1 

226 

313 

368 

1 

175 

270 

222 

8 117 92 14 


TABLE SI: Topological signature of the complex selected by simplex rate thresholding. Each cell contains a list 
of four Betti numbers, {ho, h\, h^). Eor the low rate, /o- = 0.04 Hz, To{0) encodes the correct spatial connectedness 

of the environment {bo = 1). Eor intermediate rates, 0.05 < fcr < 0.07 Hz, To{0) may occasionally break into two 
pieces {bo = 2) and for fg- > 0.10 Hz and higher, To{0) fragments into multiple components. In higher dimensions 
D > 1, To{0) contains over a hundred noncontractible topological loops. 


A. Original B. Corrected 



TABLE S2: Topological signature of the complex selected by link-rate thresholding. (A) Eor the low rate fcr = 
0.05 Hz, To{0) occasionally produces the correct topological signature {bo = b\ = bn>\ = 0, shown in bold). 
Eor intermediate rates fg- ~ 0.07 Hz, To{0) may occasionally break into two pieces {bo = 2) and produce extra 
noncontractible loops in ID. Eor high thresholds, fg- > 0.10 Hz, To{0) fragments into multiple components. However, 
the connectivity of To{0) in higher dimensions, D > 2, is correct for all cases, which implies that To{0) contracts into 
2D. (B) After applying the correction algorithms, the selected complexes acquire correct topological signature for 
fcr < 0.07 for all maps. The corresponding learning times Tmin are listed in Suppl. Table 4. 
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A. Original B. Corrected 
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no = 5 

«0 = 7 

no = 12 


no = 2 
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no = 7 
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1 31 00 

1 11 00 

1410 


1 2 

1130 
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11 
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1 2 
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1110 
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1 4 
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1 3 
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1 3 
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43 0 

1 33 0 0 

1 16 4 0 
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1 3 

1200 

1150 

1200 


TABLE S3: Topological signatures of the complexes selected by the neighbor-selection algorithm. (A) If but one 

pair of closest vertexes is selected hq = 2, 7o(^o) breaks into multiple components. For > 5, 7o(^o) has only one 
component, but path connectivity is compromised (Zpq = k ^ !)• In higher dimensions, 7“o(^o) is contractible, 
bn>i = 0- (B) After applying the correction algorithms, the selected complexes for almost all maps acquire correct 
topological signature in ID and 2D (shown in boldface) for no >1. The corresponding learning times, Lmm, are listed 
in Suppl. Table 4. 



1 

2 

3 

4 

5 

6 

7 

8 

9 

10 

mean/std 

To 

4.4 

2.7 

2.3 

2.7 

2.8 

2.7 

2.1 

3.8 

10.7 

2.7 

3.7/2.5 

me = .05) 

3.8 

2.7 

1.9 

2.7 

3.6 

3.8 

3.8 

3.8 

2.5 

3.7 

3.2/0.7 

me = .07) 

2.1 

3.8 

3.8 

2.8 

2.8 

2.3 

2.4 

3.7 

2.5 

3.8 

3.0/0.7 

Toino = 7) 

1.9 

4.6 

2.4 

2.7 

2.8 

3.8 

oo 

3.8 

oo 

3.8 

3.2/0.9 

To{no =12) 

3.8 

2.7 

1.9 

2.7 

oo 

oo 

oo 

3.8 

3.8 

oo 

3.1/0.8 


TABLE S4: The learning times, Tmin (in minutes) computed for the selected simplicial complexes To(0) and 7o(^o) 
are similar to the learning times computed via the full temporal nerve complex T. Thus, the information about 
the topological structure of the environment emerges from the cell assembly activity as fast as from the entire pool 
of place cell coactivities. However, note that the complex 7o(^o) sometimes fails to produce a finite leaning time; 
non-convergent cases are marked by oo. 


IX. SUPPLEMENTARY MOVIES ILLUSTRATING THE FIRST OF THE TEN TESTED MAPS 

Suppl. Movie 1. The grey dots represent centers of the place fields, viewed from above, similar to 
Figure 3B. The centers of the place fields that correspond to the coactive place cells are shown in red. The 
resulting activity packet moves in the environment following the simulated rat’s trajectory. 

Suppl. Movie 2. Selection of the cell assemblies by Method 1(6 - 100) and assigning readout neurons 
to the cell assemblies. 

Suppl. Movie 3. A side projection view of the activity packet propagating in the cell assembly network, 
selected via Method 1. To emphasize that the coactive place cell combinations comprise a cell assembly 
complex To(6), the corresponding place field centers are schematically connected to the readout neurons. 

Suppl. Movie 4. The same system shown in the same projection as the Figure 3B and Suppl. Movie 1. 





